#Load libraries
#General
library(deconvSeq)
library(data.table)
library(plyr)
library(dplyr)
library(tidyverse)
library(ggplot2)
library(DESeq2)
library(matrixStats)
library(ComplexUpset)
library("RColorBrewer")
library(treemapify)
library(ggtree)
library(cowplot)
library(pheatmap)
library(MDSeq)
library(matrixStats)
#Annotation
organism="org.Mm.eg.db"
library(organism, character.only = TRUE)
library("genomation")
library("GenomicRanges")
library("biomaRt")
#GSEA and ORA
library(clusterProfiler)
library(enrichplot)
library(ggnewscale)
library(rrvgo)
#Get gene mappings and biotypes
httr::set_config(httr::config(ssl_verifypeer = FALSE))
mart <- useMart('ENSEMBL_MART_ENSEMBL')
mart <- useDataset('mmusculus_gene_ensembl', mart)
annotLookup <- getBM(
mart = mart,
attributes = c(
'ensembl_gene_id',
'ensembl_transcript_id',
'external_gene_name',
'gene_biotype',
'entrezgene_id'),
uniqueRows = TRUE)
#Read in transcript and CpG island annotations
gene.obj=readTranscriptFeatures("/home/jupyter/MOUNT/references/mm10.GRCm38.96.bed", up.flank = 2000,down.flank = 2000, remove.unusual=TRUE, unique.prom = TRUE)
cpg.obj=readFeatureFlank("/home/jupyter/MOUNT/references/cpgi.mm10.bed",feature.flank.name=c("CpGi","shores"), remove.unusual = TRUE,)
bedfile<-fread("/home/jupyter/MOUNT/references/mm10.GRCm38.96.bed",header=F)
bedfile<-subset(bedfile,select=c("V1","V4","V6","V2","V3"))
colnames(bedfile)<-c("chr","feature.name","transStrand","transStart","transStop")
##Define functions
#Function for annotating methylation matrix
annotateMeth <- function(data){
gr <- makeGRangesFromDataFrame(data, keep.extra.columns=FALSE,
start.field="start",
end.field="end")
cpgAnn=annotateWithFeatureFlank(gr,
cpg.obj$CpGi,cpg.obj$shores,feature.name="CpGi",flank.name="shores")
cpgMat<-as.data.frame(getMembers(cpgAnn))
cpgMat$CpG<-ifelse(cpgMat$CpGi==1,"Island",ifelse(cpgMat$shores==1,"Shore","Other"))
data<-cbind(data,cpgMat)
#Gene annotation
geneAnn=annotateWithGeneParts(gr,gene.obj)
geneMat<-cbind(getAssociationWithTSS(geneAnn), as.data.frame(getMembers(geneAnn)))
data<-cbind(data,geneMat)
data<-merge(data,bedfile,by=c("chr","feature.name"))
#Workaround for a bug in genomation (https://github.com/BIMSBbioinfo/genomation/issues/203)
data$intron<-ifelse((data$start<=data$transStart & data$transStrand=="+"),0,
ifelse((data$start>=data$transStop & data$transStrand=="-"),0,data$intron))
data$Region<-ifelse(data$prom==1,"Promoter",ifelse(data$exon==1 | data$intron==1,"Gene body","Intergenic"))
#Filter on genomic context and significance
data<-subset(data,Region=="Promoter")
data<-merge(data,annotLookup,by.x=c("feature.name"),by.y=c("ensembl_transcript_id"))
data
}
#Function for upset plot
createUpsetPlot <- function(data) {
df<-reshape2::dcast(data, gene~Type, length)
Groups = colnames(df)[2:3]
plotUpset<-upset(df, Groups, name='Groups', width_ratio=0.1, min_size=0, min_degree=1, set_sizes=FALSE,
queries=list(
upset_query(set='Control', fill='gray'),
upset_query(set='IR', fill='magenta')
),
themes=upset_default_themes(text=element_text(size=24,face="bold"),
axis.title=element_text(face="bold",size=26),
panel.border = element_rect(colour = "black", fill=NA, size=1)),
sort_intersections_by='cardinality',
base_annotations=list('Size'=(intersection_size(counts=TRUE))),
matrix=intersection_matrix(geom=geom_point(shape='square filled',size=5)))
plotUpset
}
#Sample manifest
manifest<-fread("/home/jupyter/MOUNT/integrated/manifest.txt",header=T,sep="\t")
manifest<- manifest %>% dplyr::filter(Brain_rna!="Mmus_C57_6J_BRN_HLU_IRC_4mon_Rep5_M93") %>% dplyr::filter(Tissue=="Both" & Group %in% c("IR","Control"))
manifest[1:2,]
| Subject | Rep | Tissue | Sample | Timepoint | Group | Brain_meth | Retina_meth | Brain_rna | Retina_rna | Count | Brain_age | Retina_age |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <int> | <dbl> | <dbl> |
| M83 | Rep1 | Both | HLLC_IRC_4mon_Rep1_M83 | 4m | Control | CFG1778 | CFG1837 | Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep1_M83 | Mmus_C57_6J_RTN_HLLC_IRC_4mon_Rep1_M83 | 4 | 13.3 | 3.4 |
| M85 | Rep2 | Both | HLLC_IRC_4mon_Rep2_M85 | 4m | Control | CFG1780 | CFG1839 | Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep2_M85 | Mmus_C57_6J_RTN_HLLC_IRC_4mon_Rep2_M85 | 4 | 14.1 | NA |
#Process and merge brain methylation data
samples_202=manifest$Brain_meth
files_202=paste("GLDS-202_Epigenomics_R1_",samples_202,
"_trimmed.fq_trimmed_bismark_bt2.bismark.cov.gz",sep="")
setwd("/home/jupyter/MOUNT/integrated/bismark/glds202")
methmat_202_all = getmethmat(filnames=files_202, sample.id=samples_202, filtype="bismarkCoverage")
methmat_202_all<-annotateMeth(methmat_202_all)
#Process and merge retina methylation data
samples_203=manifest$Retina_meth
files_203=paste("GLDS-203_Epigenomics_R1_",samples_203,
"_trimmed.fq_trimmed_bismark_bt2.bismark.cov.gz",sep="")
setwd("/home/jupyter/MOUNT/integrated/bismark/glds203")
methmat_203_all = getmethmat(filnames=files_203, sample.id=samples_203, filtype="bismarkCoverage")
methmat_203_all<-annotateMeth(methmat_203_all)
manifest
| Subject | Rep | Tissue | Sample | Timepoint | Group | Brain_meth | Retina_meth | Brain_rna | Retina_rna | Count | Brain_age | Retina_age |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <int> | <dbl> | <dbl> |
| M83 | Rep1 | Both | HLLC_IRC_4mon_Rep1_M83 | 4m | Control | CFG1778 | CFG1837 | Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep1_M83 | Mmus_C57_6J_RTN_HLLC_IRC_4mon_Rep1_M83 | 4 | 13.3 | 3.4 |
| M85 | Rep2 | Both | HLLC_IRC_4mon_Rep2_M85 | 4m | Control | CFG1780 | CFG1839 | Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep2_M85 | Mmus_C57_6J_RTN_HLLC_IRC_4mon_Rep2_M85 | 4 | 14.1 | NA |
| M75 | Rep1 | Both | HLLC_IR_4mon_Rep1_M75 | 4m | IR | CFG1771 | CFG1831 | Mmus_C57_6J_BRN_HLLC_IR_4mon_Rep1_M75 | Mmus_C57_6J_RTN_HLLC_IR_4mon_Rep1_M75 | 4 | 15.3 | 2.3 |
| M80 | Rep2 | Both | HLLC_IR_4mon_Rep2_M80 | 4m | IR | CFG1775 | CFG1835 | Mmus_C57_6J_BRN_HLLC_IR_4mon_Rep2_M80 | Mmus_C57_6J_RTN_HLLC_IR_4mon_Rep2_M80 | 4 | NA | 3.7 |
| M82 | Rep4 | Both | HLLC_IR_4mon_Rep4_M82 | 4m | IR | CFG1777 | CFG1836 | Mmus_C57_6J_BRN_HLLC_IR_4mon_Rep4_M82 | Mmus_C57_6J_RTN_HLLC_IR_4mon_Rep3_M82 | 4 | 14.3 | 1.4 |
| M84 | Rep5 | Both | HLLC_IR_4mon_Rep5_M84 | 4m | IR | CFG1779 | CFG1838 | Mmus_C57_6J_BRN_HLLC_IR_4mon_Rep5_M84 | Mmus_C57_6J_RTN_HLLC_IR_4mon_Rep4_M84 | 4 | NA | 2.4 |
| M89 | Rep6 | Both | HLLC_IR_4mon_Rep6_M89 | 4m | IR | CFG1783 | CFG1842 | Mmus_C57_6J_BRN_HLLC_IR_4mon_Rep6_M89 | Mmus_C57_6J_RTN_HLLC_IR_4mon_Rep5_M89 | 4 | 11.9 | 2.8 |
| M95 | Rep3 | Both | HLLC_IRC_4mon_Rep3_M95 | 4m | Control | CFG1788 | CFG1846 | Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep3_M95 | Mmus_C57_6J_RTN_HLLC_IRC_4mon_Rep3_M95 | 4 | 12.8 | 2.2 |
| M97 | Rep4 | Both | HLLC_IRC_4mon_Rep4_M97 | 4m | Control | CFG1789 | CFG1847 | Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep4_M97 | Mmus_C57_6J_RTN_HLLC_IRC_4mon_Rep4_M97 | 4 | 14.2 | 2.5 |
#Only keep sites with >=10 reads in all samples
selectedCols <- methmat_202_all[grep("coverage", names(methmat_202_all), fixed = T)]
methmat_202_all <- methmat_202_all[apply(selectedCols, 1,
function(x) all(!is.na(x) & x>=10)), ]
selectedCols <- methmat_203_all[grep("coverage", names(methmat_203_all), fixed = T)]
methmat_203_all <- methmat_203_all[apply(selectedCols, 1,
function(x) all(!is.na(x) & x>=10)), ]
colMeans(methmat_202_all[grep("coverage", names(methmat_202_all), fixed = T)])
colMeans(methmat_203_all[grep("coverage", names(methmat_203_all), fixed = T)])
dim(methmat_202_all)
#Calculate methylation % for all sites/samples
methmat_202<-methmat_202_all
methmat_203<-methmat_203_all
rm(methmat_202_all)
rm(methmat_203_all)
for (i in 1:9){
numCs <- rlang::sym(paste("numCs", i, sep = ""))
coverage <- rlang::sym(paste("coverage", i, sep = ""))
methmat_202<-methmat_202 %>% dplyr::mutate(!!paste0("ratio_",manifest$Brain_meth[i]) := !!numCs/!!coverage)
methmat_203<-methmat_203 %>% dplyr::mutate(!!paste0("ratio_",manifest$Retina_meth[i]) := !!numCs/!!coverage)
}
#Calculate ratios (per site) and filter sites with sd<0.02
methmat_202 <- methmat_202 %>% mutate(mean = select(., starts_with("ratio_")) %>% rowMeans(na.rm = TRUE),
sd = rowSds(as.matrix(select(.,starts_with("ratio_"))),na.rm = TRUE)) %>%
filter(sd > 0.02)
methmat_203 <- methmat_203 %>% mutate(mean = select(., starts_with("ratio_")) %>% rowMeans(na.rm = TRUE),
sd = rowSds(as.matrix(select(.,starts_with("ratio_"))),na.rm = TRUE)) %>%
filter(sd > 0.02)
dim(methmat_202)
dim(methmat_203)
for (i in 1:9){
numCs <- rlang::sym(paste("numCs", i, sep = ""))
coverage <- rlang::sym(paste("coverage", i, sep = ""))
methmat_202<-methmat_202 %>% dplyr::mutate(!!coverage := ifelse(!!sym(numCs)==0,NA,!!sym(coverage)), !!numCs := ifelse(!!sym(numCs)==0,NA,!!sym(numCs)))
methmat_203<-methmat_203 %>% dplyr::mutate(!!coverage := ifelse(!!sym(numCs)==0,NA,!!sym(coverage)), !!numCs := ifelse(!!sym(numCs)==0,NA,!!sym(numCs)))
}
#Calculate weighted methylation for each gene
methmat_202_counts<- methmat_202 %>% dplyr::select(c("external_gene_name",contains(c("numCs","coverage")))) %>%
dplyr::group_by(external_gene_name) %>%
summarise(across(starts_with(c('numCs', 'coverage')), sum, na.rm=TRUE))
methmat_203_counts<- methmat_203 %>% dplyr::select(c("external_gene_name",contains(c("numCs","coverage")))) %>%
dplyr::group_by(external_gene_name) %>%
summarise(across(starts_with(c('numCs', 'coverage')), sum, na.rm=TRUE))
for (i in 1:9){
coverage <- rlang::sym(paste("coverage", i, sep = ""))
numCs <- rlang::sym(paste("numCs", i, sep = ""))
methmat_202_counts<-methmat_202_counts %>% mutate(!!paste("ratio",manifest$Brain_meth[i],sep="_") := !!numCs/!!coverage)
methmat_203_counts<-methmat_203_counts %>% mutate(!!paste("ratio",manifest$Retina_meth[i],sep="_") := !!numCs/!!coverage)
}
rm(methmat_202)
rm(methmat_203)
#Methylation ratio distribution
df_202 <- as.data.frame(gather(as.data.frame(methmat_202_counts %>% select(contains("ratio_")))))
df_202$key<-gsub("ratio_","",df_202$key)
df_203 <- as.data.frame(gather(as.data.frame(methmat_203_counts %>% select(contains("ratio_")))))
df_203$key<-gsub("ratio_","",df_203$key)
df_202<-merge(df_202,manifest,by.x=c("key"),by.y=c("Brain_meth"))
df_203<-merge(df_203,manifest,by.x=c("key"),by.y=c("Retina_meth"))
ggplot(df_202) + geom_density(aes(x=value,color=Sample))
ggplot(df_203) + geom_density(aes(x=value,color=Sample))
Warning message: “Removed 21351 rows containing non-finite values (stat_density).” Warning message: “Removed 19092 rows containing non-finite values (stat_density).”
#Compute PCs
tmp<-methmat_203_counts[,grep("ratio", names(methmat_203_counts))]
tmp[is.na(tmp)] <- 0
dfWide.pca <- prcomp(t(tmp), center = TRUE, scale. = TRUE)
dfWide.pca <- as.data.frame(dfWide.pca$x)
#PCA plots
#dfWide.pca$Tissue=rep(c("BRN","RTN"),each=length(manifest$Sample))
dfWide.pca$Group=rep(manifest$Group,times=1)
ggplot(
dfWide.pca,
aes(x = PC3, y = PC4,color=Group)
) +
scale_shape_manual(values=c(17,16)) +
geom_point(size=6, alpha = 0.7)+
geom_hline(yintercept = 0, colour = "gray65") +
geom_vline(xintercept = 0, colour = "gray65") + theme_classic()
genesCoding<-unique(as.data.frame(subset(annotLookup, gene_biotype=="protein_coding" | grepl("RNA",gene_biotype), select=c("ensembl_gene_id"))))
colnames(genesCoding) <-c("gene")
#BRN dataset (GLDS-202)
counts202<-read.delim("/home/jupyter/MOUNT/integrated/GLDS-202_rna_seq_Unnormalized_Counts.csv",sep=",",row.names=1,header=TRUE)
counts202<-counts202 %>% dplyr::filter(!grepl("ERCC",row.names(counts202))) %>%
dplyr::select(manifest$Brain_rna)
genesCounts<-as.data.frame(row.names(counts202))
colnames(genesCounts) <- c("gene")
merged<-merge(genesCounts,genesCoding,by=c("gene"))
counts202 <- counts202[merged$gene,]
coldata202<-as.data.frame(colnames(counts202))
colnames(coldata202)<-c("SampleName")
coldata202$Group<-manifest$Group
dds202 <- DESeqDataSetFromMatrix(countData = round(counts202),colData = coldata202,design =~ 1 )
ddsNofilt202 <- estimateSizeFactors(dds202)
keep <- rowSums(counts(ddsNofilt202, normalized=TRUE) >= 5 ) >= 4 & rowSums(counts(ddsNofilt202, normalized=TRUE) == 0 ) <= 4
dds202 <- ddsNofilt202[keep,]
counts202<-log2(as.data.frame(counts(dds202, normalized=TRUE))+1)
counts202$ensembl_gene_id<-rownames(counts202)
counts202<-merge(counts202,unique(subset(annotLookup,select=c("ensembl_gene_id","external_gene_name"))),by=c("ensembl_gene_id"))
#RTN dataset (GLDS-203)
counts203<-read.delim("/home/jupyter/MOUNT/integrated/GLDS-203_rna_seq_Unnormalized_Counts.csv",sep=",",row.names=1,header=TRUE)
counts203<-counts203 %>% filter(!grepl("ERCC",row.names(counts203))) %>%
dplyr::select(manifest$Retina_rna)
genesCounts<-as.data.frame(row.names(counts203))
colnames(genesCounts) <-c("gene")
merged<-merge(genesCounts,genesCoding,by=c("gene"))
counts203 <- counts203[merged$gene,]
coldata203<-as.data.frame(colnames(counts203))
colnames(coldata203)<-c("SampleName")
coldata203$Group<-manifest$Group
dds203 <- DESeqDataSetFromMatrix(countData = round(counts203),colData = coldata203,design =~ 1 )
ddsNofilt203 <- estimateSizeFactors(dds203)
keep <- rowSums(counts(ddsNofilt203, normalized=TRUE) >= 5 ) >= 4 & rowSums(counts(ddsNofilt203, normalized=TRUE) == 0 ) <= 4
dds203 <- ddsNofilt203[keep,]
counts203<-log2(as.data.frame(counts(dds203, normalized=TRUE))+1)
counts203$ensembl_gene_id<-rownames(counts203)
counts203<-merge(counts203,unique(subset(annotLookup,select=c("ensembl_gene_id","external_gene_name"))),by=c("ensembl_gene_id"))
dim(counts202)
dim(counts203)
converting counts to integer mode converting counts to integer mode
Within each exposure group, calculate Pearson correlation for each gene promoter (CpGi/s) using brain-retina matched pairs
#Merge brain and retina data on common genes
mergedMeth<-merge(methmat_202_counts %>% dplyr::select(c("external_gene_name",contains("ratio"))),
methmat_203_counts %>% dplyr::select(c("external_gene_name",contains("ratio"))),
by=c("external_gene_name"))
mergedMeth[is.na(mergedMeth)] <- 0
#Calculate correlation between brain and retina for each gene within each group
for (group in c("IR","Control")){
tmp202<-as.matrix(mergedMeth[,paste0("ratio_",subset(manifest,Group==group)$Brain_meth)])
tmp203<-as.matrix(mergedMeth[,paste0("ratio_",subset(manifest,Group==group)$Retina_meth)])
corMatTmp<-suppressWarnings(sapply(seq.int(dim(tmp202)[1]),
function(i) cor.test(tmp202[i,], tmp203[i,])))
corMatTmp<-as.data.frame(t(corMatTmp)) %>% dplyr::select(statistic,p.value,estimate,conf.int) %>%
rename_with( ~ paste0(.x,".",group))
assign(paste0("corMat_",group),corMatTmp)
}
corMatMeth<-cbind(corMat_IR,corMat_Control)
corMatMeth$gene<-mergedMeth$external_gene_name
posCorMeth[1,]
| gene | Type | |
|---|---|---|
| <chr> | <chr> | |
| 1 | 1110006O24Rik | IR |
#Positively correlated
posCorMeth<-unique(rbind(as.data.frame(corMatMeth %>% dplyr::filter(p.value.IR<0.05 & estimate.IR>0)
%>% dplyr::select(gene) %>% mutate(Type="IR")),
as.data.frame(corMatMeth %>% dplyr::filter(p.value.Control<0.05 & estimate.Control>0)
%>% dplyr::select(gene) %>% mutate(Type="Control"))))
#Negatively correlated
negCorMeth<-unique(rbind(as.data.frame(corMatMeth %>% dplyr::filter(p.value.IR<0.05 & estimate.IR<0)
%>% dplyr::select(gene) %>% mutate(Type="IR")),
as.data.frame(corMatMeth %>% dplyr::filter(p.value.Control<0.05 & estimate.Control<0)
%>% dplyr::select(gene) %>% mutate(Type="Control"))))
options(repr.plot.width=6, repr.plot.height=8)
plot_grid(createUpsetPlot(posCorMeth), createUpsetPlot(negCorMeth),ncol=1)
posCorMeth$assay<-"RRBS"
negCorMeth$assay<-"RRBS"
Using Type as value column: use value.var to override. Using Type as value column: use value.var to override.
mergedRna<-merge(counts202,counts203,by="ensembl_gene_id")
for (group in c("IR","Control")){
tmp202<-as.matrix(t(apply(mergedRna[,subset(manifest,Group==group)$Brain_rna],1,scale)))
tmp203<-as.matrix(t(apply(mergedRna[,subset(manifest,Group==group)$Retina_rna],1,scale)))
corMatTmp<-suppressWarnings(sapply(seq.int(dim(tmp202)[1]),
function(i) cor.test(tmp202[i,], tmp203[i,])))
corMatTmp<-as.data.frame(t(corMatTmp)) %>% dplyr::select(statistic,p.value,estimate,conf.int) %>%
rename_with( ~ paste0(.x,".",group))
assign(paste0("corMatRna_",group),corMatTmp)
}
corMatRna<-cbind(corMatRna_IR,corMatRna_Control)
corMatRna$ensembl_gene_id<-mergedRna$ensembl_gene_id
corMatRna<-merge(corMatRna,unique(subset(annotLookup,select=c("ensembl_gene_id","external_gene_name"))),by=c("ensembl_gene_id"))
corMatRna <- corMatRna %>% rename(gene=external_gene_name)
#Positively correlated
posCorRna<-unique(rbind(as.data.frame(corMatRna %>% dplyr::filter(p.value.IR<0.05 & estimate.IR>0)
%>% dplyr::select(gene) %>% mutate(Type="IR")),
as.data.frame(corMatRna %>% dplyr::filter(p.value.Control<0.05 & estimate.Control>0)
%>% dplyr::select(gene) %>% mutate(Type="Control"))))
#Negatively correlated
negCorRna<-unique(rbind(as.data.frame(corMatRna %>% dplyr::filter(p.value.IR<0.05 & estimate.IR<0)
%>% dplyr::select(gene) %>% mutate(Type="IR")),
as.data.frame(corMatRna %>% dplyr::filter(p.value.Control<0.05 & estimate.Control<0)
%>% dplyr::select(gene) %>% mutate(Type="Control"))))
options(repr.plot.width=6, repr.plot.height=8)
plot_grid(createUpsetPlot(posCorRna), createUpsetPlot(negCorRna),nrow=2)
posCorRna$assay<-"RNA"
negCorRna$assay<-"RNA"
Using Type as value column: use value.var to override. Using Type as value column: use value.var to override.
allPosCor<-rbind(posCorMeth,posCorRna)
allPosCor$Type<-paste(allPosCor$assay,allPosCor$Type,sep="_")
df<-reshape2::dcast(allPosCor, gene~Type, length)
df[1,]
Groups = colnames(df)[2:5]
plotUpset<-upset(df, Groups, name='Groups', width_ratio=0.1, min_size=0, min_degree=1, set_sizes=FALSE,
queries=list(
upset_query(set='RNA_Control', fill='lightblue'),
upset_query(set='RNA_IR', fill='orange'),
upset_query(set='RRBS_Control', fill='blue'),
upset_query(set='RRBS_IR', fill='brown')
),
themes=upset_default_themes(text=element_text(size=24,face="bold"),
axis.title=element_text(face="bold",size=26),
panel.border = element_rect(colour = "black", fill=NA, size=1)),
sort_intersections_by='cardinality',
base_annotations=list('Size'=(intersection_size(counts=TRUE))),
matrix=intersection_matrix(geom=geom_point(shape='square filled',size=5)))
options(repr.plot.width=15, repr.plot.height=7)
plotUpset
dfSum = df %>% mutate(sum = select(., !starts_with("gene") & !contains("Control")) %>% rowSums(na.rm = TRUE))
Using assay as value column: use value.var to override.
| gene | RNA_Control | RNA_IR | RRBS_Control | RRBS_IR | |
|---|---|---|---|---|---|
| <chr> | <int> | <int> | <int> | <int> | |
| 1 | 0610010K14Rik | 0 | 0 | 1 | 0 |
subset(dfSum, sum>0 & RRBS_Control==0 & RNA_Control==0 & RNA_IR==1 & RRBS_IR==1)$gene
Within each exposure group for each tissue type, calculate Pearson correlation for each gene promoter (CpGi/s) using RRBS-RNA-seq matched pairs
methmat_203_counts[is.na(methmat_203_counts)] <- 0
methmat_202_counts[is.na(methmat_202_counts)] <- 0
merged203<-merge(counts203,methmat_203_counts,by="external_gene_name")
for (group in c("IR","Control")){
tmpMeth<-as.matrix(merged203[,paste0("ratio_",subset(manifest,Group==group)$Retina_meth)])
tmpRna<-as.matrix(t(apply(merged203[,subset(manifest,Group==group)$Retina_rna],1,scale)))
corMat<-suppressWarnings(sapply(seq.int(dim(tmpMeth)[1]),
function(i) cor.test(tmpMeth[i,], tmpRna[i,])))
corMat<-as.data.frame(t(corMat)) %>% dplyr::select(statistic,p.value,estimate,conf.int) %>%
rename_with( ~ paste0(.x,".",group))
assign(paste0("corMatRetina_",group),corMat)
}
merged202<-merge(counts202,methmat_202_counts,by="external_gene_name")
for (group in c("IR","Control")){
samples_rna<-subset(manifest,Group==group)$Brain_rna
samples_meth<-subset(manifest,Group==group)$Brain_meth
tmpMeth<-as.matrix(merged202[,paste0("ratio_",samples_meth)])
tmpRna<-as.matrix(t(apply(merged202[,samples_rna],1,scale)))
corMat<-suppressWarnings(sapply(seq.int(dim(tmpMeth)[1]),
function(i) cor.test(tmpMeth[i,], tmpRna[i,])))
corMat<-as.data.frame(t(corMat)) %>% dplyr::select(statistic,p.value,estimate,conf.int) %>%
rename_with( ~ paste0(.x,".",group))
assign(paste0("corMatBrain_",group),corMat)
}
samples_rna<-subset(manifest,Group==group)$Brain_rna
samples_meth<-subset(manifest,Group==group)$Brain_meth
ggplot(subset(merged202, !!subset(manifest,Group==group)$Brain_rna[2] > 5),aes(x=get(subset(manifest,Group==group)$Brain_rna[2]),
y=get(paste0("ratio_",samples_meth[2])))) + geom_point() + geom_smooth(method="lm")
`geom_smooth()` using formula 'y ~ x'
corMatBrain<-cbind(corMatBrain_IR,corMatBrain_Control)
corMatRetina<-cbind(corMatRetina_IR,corMatRetina_Control)
corMatBrain$gene<-merged202$external_gene_name
corMatRetina$gene<-merged203$external_gene_name
dim(corMatBrain)
dim(corMatRetina)
t(subset(merged202, external_gene_name=="Fktn"))
| 3864 | |
|---|---|
| external_gene_name | Fktn |
| ensembl_gene_id | ENSMUSG00000028414 |
| Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep1_M83 | 9.753355 |
| Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep2_M85 | 9.834692 |
| Mmus_C57_6J_BRN_HLLC_IR_4mon_Rep1_M75 | 9.672554 |
| Mmus_C57_6J_BRN_HLLC_IR_4mon_Rep2_M80 | 9.708237 |
| Mmus_C57_6J_BRN_HLLC_IR_4mon_Rep4_M82 | 9.369445 |
| Mmus_C57_6J_BRN_HLLC_IR_4mon_Rep5_M84 | 9.746524 |
| Mmus_C57_6J_BRN_HLLC_IR_4mon_Rep6_M89 | 9.762727 |
| Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep3_M95 | 10.03492 |
| Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep4_M97 | 9.790137 |
| numCs1 | 2 |
| numCs2 | 3 |
| numCs3 | 1 |
| numCs4 | 74 |
| numCs5 | 67 |
| numCs6 | 118 |
| numCs7 | 5 |
| numCs8 | 5 |
| numCs9 | 1 |
| coverage1 | 106 |
| coverage2 | 115 |
| coverage3 | 59 |
| coverage4 | 1350 |
| coverage5 | 555 |
| coverage6 | 2688 |
| coverage7 | 220 |
| coverage8 | 240 |
| coverage9 | 40 |
| ratio_CFG1778 | 0.01886792 |
| ratio_CFG1780 | 0.02608696 |
| ratio_CFG1771 | 0.01694915 |
| ratio_CFG1775 | 0.05481481 |
| ratio_CFG1777 | 0.1207207 |
| ratio_CFG1779 | 0.04389881 |
| ratio_CFG1783 | 0.02272727 |
| ratio_CFG1788 | 0.02083333 |
| ratio_CFG1789 | 0.025 |
#Positively correlated between RRBS and RNA-Seq
posCorRetina<-unique(rbind(as.data.frame(corMatRetina %>% dplyr::filter(p.value.IR<0.05 & estimate.IR>0)
%>% dplyr::select(gene) %>% mutate(Type="IR")),
as.data.frame(corMatRetina %>% dplyr::filter(p.value.Control<0.05 & estimate.Control>0)
%>% dplyr::select(gene) %>% mutate(Type="Control"))))
#Negatively correlated between RRBS and RNA-Seq
negCorRetina<-unique(rbind(as.data.frame(corMatRetina %>% dplyr::filter(p.value.IR<0.05 & estimate.IR<0)
%>% dplyr::select(gene) %>% mutate(Type="IR")),
as.data.frame(corMatRetina %>% dplyr::filter(p.value.Control<0.05 & estimate.Control<0)
%>% dplyr::select(gene) %>% mutate(Type="Control"))))
options(repr.plot.width=6, repr.plot.height=8)
plot_grid(createUpsetPlot(posCorRetina), createUpsetPlot(negCorRetina),ncol=1)
posCorRetina$tissue<-"RTN"
negCorRetina$tissue<-"RTN"
Using Type as value column: use value.var to override. Using Type as value column: use value.var to override.
#Positively correlated between RRBS and RNA-Seq
posCorBrain<-unique(rbind(as.data.frame(corMatBrain %>% dplyr::filter(p.value.IR<0.05 & estimate.IR>0)
%>% dplyr::select(gene) %>% mutate(Type="IR")),
as.data.frame(corMatBrain %>% dplyr::filter(p.value.Control<0.05 & estimate.Control>0)
%>% dplyr::select(gene) %>% mutate(Type="Control"))))
#Negatively correlated between RRBS and RNA-Seq
negCorBrain<-unique(rbind(as.data.frame(corMatBrain %>% dplyr::filter(p.value.IR<0.05 & estimate.IR<0)
%>% dplyr::select(gene) %>% mutate(Type="IR")),
as.data.frame(corMatBrain %>% dplyr::filter(p.value.Control<0.05 & estimate.Control<0)
%>% dplyr::select(gene) %>% mutate(Type="Control"))))
options(repr.plot.width=6, repr.plot.height=8)
plot_grid(createUpsetPlot(posCorBrain), createUpsetPlot(negCorBrain),ncol=1)
posCorBrain$tissue<-"BRN"
negCorBrain$tissue<-"BRN"
Using Type as value column: use value.var to override. Using Type as value column: use value.var to override.
allPosCor<-rbind(negCorBrain,negCorRetina)
allPosCor$Type<-paste(allPosCor$tissue,allPosCor$Type,sep="_")
df<-reshape2::dcast(allPosCor, gene~Type, length)
Groups = colnames(df)[2:5]
plotUpset2<-upset(df, Groups, name='Groups', width_ratio=0.1, min_size=0, min_degree=1, set_sizes=FALSE,
queries=list(
upset_query(set='BRN_Control', fill='lightblue'),
upset_query(set='BRN_IR', fill='orange'),
upset_query(set='RTN_Control', fill='blue'),
upset_query(set='RTN_IR', fill='brown')
),
themes=upset_default_themes(text=element_text(size=24,face="bold"),
axis.title=element_text(face="bold",size=26),
panel.border = element_rect(colour = "black", fill=NA, size=1)),
sort_intersections_by='cardinality',
base_annotations=list('Size'=(intersection_size(counts=TRUE))),
matrix=intersection_matrix(geom=geom_point(shape='square filled',size=5)))
options(repr.plot.width=15, repr.plot.height=7)
plotUpset2
dfSum = df %>% mutate(sum = select(., !starts_with("gene") & !contains("Control")) %>% rowSums(na.rm = TRUE))
Using tissue as value column: use value.var to override.
| gene | BRN_Control | BRN_IR | RTN_Control | RTN_IR | |
|---|---|---|---|---|---|
| <chr> | <int> | <int> | <int> | <int> | |
| 1 | 0610040J01Rik | 0 | 1 | 0 | 0 |
subset(dfSum, sum>1 & BRN_Control==0 & RTN_Control==0)
| gene | BRN_Control | BRN_IR | RTN_Control | RTN_IR | sum | |
|---|---|---|---|---|---|---|
| <chr> | <int> | <int> | <int> | <int> | <dbl> | |
| 16 | 5033430I15Rik | 0 | 1 | 0 | 1 | 2 |
| 339 | Eln | 0 | 1 | 0 | 1 | 2 |
| 394 | Fgf11 | 0 | 1 | 0 | 1 | 2 |
| 401 | Fktn | 0 | 1 | 0 | 1 | 2 |
| 413 | G3bp2 | 0 | 1 | 0 | 1 | 2 |
| 454 | Gm28778 | 0 | 1 | 0 | 1 | 2 |
| 668 | Ncs1 | 0 | 1 | 0 | 1 | 2 |
| 797 | Ppm1e | 0 | 1 | 0 | 1 | 2 |
| 924 | Sgsh | 0 | 1 | 0 | 1 | 2 |
options(repr.plot.width=20, repr.plot.height=8)
plot_grid(plotUpset,plotUpset2,nrow=1,labels=c("A","B"))
setwd("/home/jupyter/glds202_203/Pathways")
cases<-c("IR","Control")
for (j in 1:length(cases)){
title<-cases[j]
tmpData<-subset(as.data.frame(corMatRetina),select=c(paste0("estimate.",title),paste0("p.value.",title),"gene"))
colnames(tmpData)<-c("estimate","pvalue","gene")
tmpData<-subset(tmpData,!is.na(estimate))
tmpData$pvalue<-ifelse(tmpData$pvalue==0, 0.0000000001, tmpData$pvalue)
#tmpData$stat<-sign(as.numeric(tmpData$estimate))* -log10(as.numeric(tmpData$pvalue))
tmpData$stat<-as.numeric(tmpData$estimate) * -log10(as.numeric(tmpData$pvalue))
tmpData<-subset(tmpData,select=c("stat","gene"))
tmpData<-subset(tmpData,!is.na(stat))
rankMetric <- as.numeric(tmpData$stat)
names(rankMetric) <- tmpData$gene
rankMetric <- sort(rankMetric, decreasing = TRUE)
####GO####
gseaGO <- gseGO(geneList=rankMetric,
ont ="BP",
keyType="SYMBOL",
pvalueCutoff = 0.25,
by = "fgsea",
pAdjustMethod = "BH",
OrgDb = org.Mm.eg.db)
if(dim(gseaGO)[1]>0){
gseaRes<-data.frame(gseaGO)
gseaRes$Dataset<-"RTN"
gseaRes$Group<-title
write.table(gseaRes,paste0("RTN",title,"_clusterProfiler.GO.BP.all.symbol.tsv"),sep="\t",quote=FALSE,row.names=F)
}
}
preparing geneSet collections... GSEA analysis... Warning message in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : “There are duplicate gene names, fgsea may produce unexpected results.” no term enriched under specific pvalueCutoff... preparing geneSet collections... GSEA analysis... Warning message in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : “There are duplicate gene names, fgsea may produce unexpected results.” no term enriched under specific pvalueCutoff...
corTissue<-rbind(posCorRetina %>% mutate(Sign="pos"), negCorRetina %>% mutate(Sign="neg"),
posCorBrain %>% mutate(Sign="pos"), negCorBrain %>% mutate(Sign="neg"))
corAssay<-rbind(posCorMeth %>% mutate(Sign="pos"), negCorMeth %>% mutate(Sign="neg"),
posCorRna %>% mutate(Sign="pos"), negCorRna %>% mutate(Sign="neg"))
corTissue[1,]
subset(ddply(unique(subset(corTissue,Type=="IR",select=c("gene","tissue"))), .(gene), nrow),V1>1)$gene
subset(ddply(unique(subset(corTissue,Type=="Control",select=c("gene","tissue"))), .(gene), nrow),V1>1)$gene
| gene | Type | tissue | Sign | |
|---|---|---|---|---|
| <chr> | <chr> | <chr> | <chr> | |
| 1 | 1700029J07Rik | IR | RTN | pos |
diffVarBrain<-fread("/home/jupyter/glds202_203/DiffVar/BRN_IR_vs_Control.tsv",header=T,sep="\t")
sort(subset(diffVarBrain,
external_gene_name %in% c(subset(negCorBrain,Type=="IR")$gene,subset(posCorBrain,Type=="IR")$gene) &
(FDR.dispersion<0.05 | FDR.mean<0.05))$external_gene_name)
tmp1<-subset(counts202,grepl("Cul1",external_gene_name))
tmp2<-subset(counts203,grepl("Cul1",external_gene_name))
options(repr.plot.width=10, repr.plot.height=8)
tmp1Long<-melt(subset(tmp1,select=-c(ensembl_gene_id))) %>%
mutate(group=ifelse(grepl("HLLC_IRC",variable),"Control",ifelse(grepl("HLLC_IR_",variable),"IR",
ifelse(grepl("HLU_IRC",variable),"HLU","Combined"))),
timepoint=ifelse(grepl("4mon",variable),"4m",ifelse(grepl("1mon",variable),"1m","7d")))
ggplot(tmp1Long) + geom_boxplot(aes(x=group,y=value,color=group)) + xlab("gene") + ylab("Normalized counts") +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + facet_wrap(.~external_gene_name,scale="free") + ggtitle("Brain")
tmp2Long<-melt(subset(tmp2,select=-c(ensembl_gene_id))) %>%
mutate(group=ifelse(grepl("HLLC_IRC",variable),"Control",ifelse(grepl("HLLC_IR_",variable),"IR",
ifelse(grepl("HLU_IRC",variable),"HLU","Combined"))),
timepoint=ifelse(grepl("4mon",variable),"4m",ifelse(grepl("1mon",variable),"1m","7d")))
ggplot(tmp2Long) + geom_boxplot(aes(x=group,y=value,color=group)) + xlab("gene") + ylab("Normalized counts") +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + facet_wrap(.~external_gene_name,scale="free") + ggtitle("Retina")
Warning message in melt(subset(tmp1, select = -c(ensembl_gene_id))): “The melt generic in data.table has been passed a data.frame and will attempt to redirect to the relevant reshape2 method; please note that reshape2 is deprecated, and this redirection is now deprecated as well. To continue using melt methods from reshape2 while both libraries are attached, e.g. melt.list, you can prepend the namespace like reshape2::melt(subset(tmp1, select = -c(ensembl_gene_id))). In the next version, this warning will become an error.” Using external_gene_name as id variables Warning message in melt(subset(tmp2, select = -c(ensembl_gene_id))): “The melt generic in data.table has been passed a data.frame and will attempt to redirect to the relevant reshape2 method; please note that reshape2 is deprecated, and this redirection is now deprecated as well. To continue using melt methods from reshape2 while both libraries are attached, e.g. melt.list, you can prepend the namespace like reshape2::melt(subset(tmp2, select = -c(ensembl_gene_id))). In the next version, this warning will become an error.” Using external_gene_name as id variables
diffVarRetina<-fread("/home/jupyter/glds202_203/DiffVar/RTN_IR_vs_Control.tsv",header=T,sep="\t")
diffVarRetina<-subset(diffVarRetina,select=c("external_gene_name"), FDR.dispersion<0.1)
for (i in 1:dim(subset(manifest))[1]){
tmp<-subset(merged203,select=c(subset(manifest)$Retina_rna[i],
paste0("ratio_",subset(manifest)$Retina_meth[i])),
external_gene_name %in% diffVarRetina$external_gene_name)
colnames(tmp)<-c("Brain","Retina")
if(i==1){
longRna<-tmp
}
else{
longRna<-rbind(longRna,tmp)
}
}
ggplot(longRna,aes(x=Brain,y=Retina)) + geom_point(size=5) + geom_smooth(method="lm")
`geom_smooth()` using formula 'y ~ x'
manifest[1,]
| Subject | Rep | Tissue | Sample | Timepoint | Group | Brain_meth | Retina_meth | Brain_rna | Retina_rna | Count | Brain_age | Retina_age |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <chr> | <int> | <dbl> | <dbl> |
| M83 | Rep1 | Both | HLLC_IRC_4mon_Rep1_M83 | 4m | Control | CFG1778 | CFG1837 | Mmus_C57_6J_BRN_HLLC_IRC_4mon_Rep1_M83 | Mmus_C57_6J_RTN_HLLC_IRC_4mon_Rep1_M83 | 4 | 13.3 | 3.4 |
#Network of genes with correlation and differential expression/variability
library("STRINGdb")
string_db <- STRINGdb$new(version="11.5", species=10090,score_threshold=400,input_directory="")
diffVarBrain<-fread("/home/jupyter/glds202_203/DiffVar/BRN_IR_vs_Control.tsv",header=T,sep="\t")
diffVarBrain<-subset(diffVarBrain,select=c("external_gene_name"), FDR.dispersion<1 | FDR.mean<1)
diffVarBrain<-unique(subset(diffVarBrain,
external_gene_name %in% c(subset(posCorBrain,Type=="Control")$gene,
subset(negCorBrain,Type=="Control")$gene,
subset(posCorBrain,Type=="IR")$gene,
subset(negCorBrain,Type=="IR")$gene),
select=c("external_gene_name")))
length(diffVarBrain$external_gene_name)
cat(paste0(sort(diffVarBrain$external_gene_name),sep=","))
cat("\n")
diffVarRetina<-fread("/home/jupyter/glds202_203/DiffVar/RTN_IR_vs_Control.tsv",header=T,sep="\t")
diffVarRetina<-subset(diffVarRetina,select=c("external_gene_name"), FDR.dispersion<1 | FDR.mean<1)
diffVarRetina<-unique(subset(diffVarRetina,
external_gene_name %in% c(subset(posCorRetina,Type=="Control")$gene,
subset(negCorRetina,Type=="Control")$gene,
subset(posCorRetina,Type=="IR")$gene,
subset(negCorRetina,Type=="IR")$gene),
select=c("external_gene_name")))
length(diffVarRetina$external_gene_name)
cat(paste0(sort(diffVarRetina$external_gene_name),sep=","))
df<-merge(diffVarBrain,diffVarRetina,by=c("external_gene_name"))
colnames(df)<-c("gene")
df<-diffVarBrain
colnames(df)<-c("gene")
diffVarBrain$group<-"BRN_IR"
diffVarRetina$group<-"RTN_IR"
df<-rbind(diffVarBrain,diffVarRetina)
0610040J01Rik, 1500009L16Rik, 1700030C10Rik, 1700034H15Rik, 1700102P08Rik, 2410018L13Rik, 2510009E07Rik, 2700049A03Rik, 2810455O05Rik, 2810459M11Rik, 3010001F23Rik, 3110009E18Rik, 4632411P08Rik, 4930447C04Rik, 4933427D14Rik, 5033430I15Rik, 5330434G04Rik, 5830408C22Rik, 6430590A07Rik, 9030624G23Rik, 9130024F11Rik, 9330111N05Rik, 9330188P03Rik, 9530068E07Rik, A330102I10Rik, A630072M18Rik, A930002C04Rik, Aagab, Abca5, Abcb1b, Abcb4, Abcb8, Abcf1, Abcg1, Abracl, Acadvl, Ache, Ackr3, Acot1, Acot13, Acot4, Acpp, Actn4, Actr1a, Actr2, Adam11, Adam23, Adamtsl5, Adar, Adgra2, Ado, Adra1d, Adra2c, Adssl1, Afg1l, Afg3l2, Afmid, Agbl2, Ago4, AI837181, Akap13, Akr7a5, Alcam, Alg3, Alkal2, Alkbh1, Alkbh5, Alox12b, Anapc15, Anapc2, Ankrd13a, Ankrd13b, Ankrd40, Ankrd45, Ap1s2, Ap3m2, Ap4e1, Ap5s1, Apbb1, Apcdd1, Aplp1, Arcn1, Arhgap28, Arhgap42, Arhgef10l, Arhgef33, Arhgef4, Arhgef5, Arih2, Arl10, Arl6ip4, Arl6ip5, Arl8b, Arrdc3, Arxes1, Arxes2, Asb8, Asnsd1, Ass1, Atad2b, Atg16l2, Atp2b2, Atp5g2, Atp5o, Atpaf1, AU040320, Avpi1, B230219D22Rik, B4galnt4, B4gat1, B630019K06Rik, Bahcc1, Baz2a, Bbs4, Bbs5, Bbs7, BC029722, BC037039, Bcap29, Bcat2, Bcdin3d, Bcl7c, Bean1, Bend4, Bend5, Bicdl1, Bin2, Birc6, Blvra, Bmpr1b, Boc, Brd1, Brd7, Bri3, Brms1l, Bsn, Btbd8, Btc, C030037D09Rik, C1ql3, C2, C230014O12Rik, C2cd2l, C330011M18Rik, Cactin, Cadps, Camk2d, Camk2n2, Caml, Camsap2, Capns1, Capza2, Car11, Car4, Caskin1, Cbr2, Ccdc113, Ccdc158, Ccdc160, Ccdc162, Ccdc22, Ccdc88c, Ccdc9, Ccdc9b, Ccl27a, Ccn1, Ccnk, Cdc14a, Cdc16, Cdc40, Cdh4, Cdh9, Cdk13, Cdk20, Cdk5rap2, Cdkn1b, Cdkn1c, Cdkn2d, Cds2, Cdv3, Celsr3, Cep63, Cep70, Cerox1, Cert1, Cfap298, Cfap36, Cfap418, Cfap57, Cfap61, Cfap97d2, Cftr, Chchd3, Chka, Chordc1, Chst10, Chst15, Chsy3, Chtop, Cilk1, Cipc, Ckap2l, Clcn4, Clcn5, Cldn1, Cldn12, Clspn, Cmtr2, Cnih1, Cnst, Cntnap2, Cntrl, Col9a2, Coq10b, Coq7, Cox17, Cox5b, Cox6c, Cpeb1, Cpne5, Cpq, Cpsf7, Crb2, Crem, Crtc1, Csf1r, Csnk1g3, Ctf1, Ctnnal1, Ctnnb1, Ctsh, Ctsz, Cul5, Cutc, Cxcl12, Cyb5a, Cycs, Cygb, Cyp46a1, Cyth3, D130020L05Rik, D430036J16Rik, Daam2, Dach1, Dcaf10, Dcakd, Dcbld1, Ddit4l, Ddx46, Ddx49, Dennd4c, Desi1, Dhtkd1, Dhx29, Diablo, Dido1, Dip2b, Dip2c, Dis3, Dlat, Dlg4, Dlg5, Dmpk, Dnaaf3, Dnah9, Dnaja3, Dnajb2, Dnajc15, Dnajc3, Dnmbp, Dnmt3a, Dok7, Dolpp1, Dpp6, Dpy19l1, Dpysl3, Dpysl5, Drd1, Dtx4, Dusp16, Dusp27, Dusp8, Dvl3, Dym, Dyrk1b, Ecpas, Edem2, Eef1g, Efcab12, Eid2b, Eif1b, Eif3e, Eif3f, Eif3l, Eif5a, Eif5a2, Eln, Emc10, Eml1, Emp3, En2, Epm2a, Epm2aip1, Eps8, Erbb3, Erich5, Errfi1, Esyt1, Ets1, Etv1, Exoc3l, Exoc3l4, Exosc3, Exosc5, Exosc8, Fads2, Faf2, Fahd1, Faim2, Fam107a, Fam110a, Fam117b, Fam118b, Fam131c, Fam136a, Fam163b, Fam171a1, Fam172a, Fam189a1, Fam207a, Fam234b, Fam53c, Fbln1, Fbn2, Fbxl17, Fbxl20, Fbxl3, Fbxo21, Fbxo6, Fbxw7, Fbxw9, Fdps, Fem1b, Fgd3, Fgd4, Fgf11, Fgf9, Fgfbp3, Fhl3, Fign, Fkbp1b, Fktn, Flcn, Fn3krp, Frk, Fsd1, Fyb, Fzd7, G3bp2, Gabrg3, Gadd45gip1, Galc, Galnt7, Galnt9, Garnl3, Gas7, Gas8, Gatad2b, Gatd1, Gatm, Gcdh, Get3, Gfer, Ggn, Ggta1, Ghsr, Gimap1, Gins4, Glra1, Glt1d1, Gltp, Gm10545, Gm11613, Gm11642, Gm12216, Gm12743, Gm13889, Gm14325, Gm15609, Gm17435, Gm20387, Gm26559, Gm26871, Gm27003, Gm28050, Gm28778, Gm29394, Gm34552, Gm35394, Gm41760, Gm4285, Gm43254, Gm45605, Gm45823, Gm4673, Gm48623, Gm50394, Gm5602, Gm5784, Gm9828, Gm9856, Gnat1, Gng2, Golgb1, Gorasp1, Gpn1, Gpr137, Gpr139, Gpr3, Gpx7, Grb7, Gria2, Grik3, Grk2, Gsn, Gtf2e2, Gtf2h4, Gyg, Gzmm, H1f5, H3c4, H3f3b, H6pd, Haus1, Haus4, Haus5, Haus7, Hcrtr1, Hdgf, Hdgfl3, Herpud2, Hhip, Hic2, Hid1, Hif1an, Higd1a, Hint2, Hmces, Hmgcr, Hnrnph3, Homez, Hook3, Hsd17b7, Hspa8, Hspb8, Ids, Idua, Igbp1, Igfbp7, Igsf3, Il10ra, Il20ra, Il33, Impact, Ino80, Inpp5a, Ints11, Iqgap1, Ireb2, Isyna1, Itga8, Itga9, Itgb4, Jazf1, Jkamp, Jtb, Kat6b, Katna1, Katnb1, Kbtbd3, Kcnh5, Kcnip1, Kcnk1, Kcnq4, Kcp, Kctd6, Kdelr2, Kdm2a, Kdm8, Keap1, Khk, Kif1b, Kif1c, Kif3a, Kifbp, Kifc3, Kirrel2, Klhl1, Klhl14, Klhl35, Kpnb1, Krt12, Kyat3, L2hgdh, Lamc1, Larp1, Lars2, Lcmt2, Lcorl, Lgals8, Lgr6, Lima1, Limd1, Lin7c, Llgl1, Lman1, Lmtk3, Lonrf3, Loxl3, Lpcat3, Lpin2, Lrfn4, Lrp1b, Lrrc28, Lrrc49, Lrrfip2, Lrrn2, Lsm10, Lzts1, Lzts2, Macir, Mak16, Malat1, Mall, Man2a1, Map2k3, Map4k2, Map7d2, Map9, Mapk3, Mapk8, Mapre2, Marchf3, Marf1, Mast1, Mcc, Mcf2l, Mcfd2, Mcrip2, Mdm1, Med10, Med15, Med23, Meg3, Mettl5, Mfsd14b, Mfsd6, Mfsd9, Mgat5, Mgmt, Mgrn1, Mical2, Mid1, Mier1, Mindy2, Mios, Mipol1, Mkrn1, Mllt6, Mms22l, Morc2a, Mospd1, Moxd1, Mpnd, Mpp2, Mpst, Mrc2, Mroh5, Mrpl3, Msl3, Msra, Mta3, Mtg2, Mtss2, Mtus2, Mus81, Mxd4, Myadml2, Mybl2, Mycl, Myh9, Myo1d, Myo1g, Myo9b, Naa10, Naa50, Naa60, Nalcn, Nanos1, Nanp, Napa, Nav1, Ncan, Ncapd2, Nck1, Ncoa2, Ncoa6, Ncor1, Ncs1, Ndst3, Ndufa1, Ndufa12, Ndufa8, Ndufb5, Ndufs6, Necap2, Neil1, Neurod2, Nexn, Nfat5, Nfatc3, Nfatc4, Nfkbie, Nfya, Nfyb, Nhlrc2, Niban2, Nim1k, Nipsnap2, Noct, Nog, Nopchap1, Nova1, Npas4, Npl, Nploc4, Nptn, Nr3c2, Nr4a2, Nrg4, Nrgn, Nrxn2, Nsd1, Nsmce1, Ntmt1, Nubp1, Nudt16l1, Nup188, Nxn, Oaf, Obscn, Olfm1, Olfml2b, Olfml3, Oprm1, Optn, Osbpl1a, Osbpl3, Osbpl6, Ostf1, Otud1, Ovgp1, Oxsr1, P4ha3, P4hb, Pabpc1, Pacs2, Pafah1b1, Pak5, Pals2, Pamr1, Pank1, Papola, Paqr8, Paqr9, Parp1, Parp14, Parp16, Parp3, Pbxip1, Pcdhga1, Pcdhga9, Pcdhgb4, Pcdhgb6, Pced1a, Pck2, Pcnx2, Pcsk6, Pdap1, Pdcd2, Pde6d, Pde7a, Pdgfb, Pdlim2, Pdp1, Pea15a, Per3, Pex11a, Pex12, Pgam2, Pgpep1, Phactr1, Phactr2, Phc1, Phc2, Phc3, Phf13, Phf23, Phpt1, Pias2, Picalm, Pigq, Pigv, Pik3cd, Pim2, Pinx1, Pip4k2b, Pitpnm1, Pkp2, Pla2g12a, Plce1, Plekha4, Plekhd1, Plekhg5, Plekhh1, Plekhm3, Plppr1, Pls1, Pm20d2, Pnma8b, Pno1, Poglut1, Poldip2, Pole, Poll, Polr3a, Polr3b, Ppa2, Ppara, Ppcs, Pphln1, Ppia, Ppic, Ppip5k1, Ppm1e, Ppm1k, Ppp1r10, Ppp1r16b, Ppp1r3f, Ppp5c, Ppt1, Prcc, Prkab2, Prkcg, Prkch, Prkdc, Prkn, Prps2, Psmb10, Psph, Psrc1, Ptbp1, Ptcd2, Pter, Ptpn1, Ptpn13, Ptpn18, Ptprm, Ptprn, Ptpro, Puf60, Pwp1, Pxdc1, Pxdn, Pxmp2, Qrich1, Qrsl1, Rab27a, Rab33a, Rab40b, Rab42, Rabggta, Rabl6, Rad9a, Radil, Ramp1, Ranbp17, Ranbp3l, Rap1gds1, Rasa4, Rasgrp1, Rassf2, Rassf8, Rbm26, Rbm4, Rcan3, Rcbtb2, Reck, Reep4, Retreg1, Rfxap, Rgcc, Rgs7, Rgs7bp, Rlf, Rmc1, Rmi1, Rnf114, Rnf144a, Rnf166, Rnf170, Rnf182, Rnf185, Rnf217, Rnft2, Robo2, Rogdi, Rpgrip1, Rpia, Rpl12, Rpl15, Rpl24, Rpl26, Rpl29, Rpl34, Rpl35a, Rpl41, Rpl7, Rpl8, Rplp0, Rplp2, Rpn1, Rps6ka3, Rrbp1, Rreb1, Rrp8, Rrs1, Rsc1a1, Rwdd3, Rxfp3, S100a11, S1pr5, Sae1, Samd11, Samd15, Sart1, Scfd2, Scly, Scn2b, Scnn1a, Sdf2, Sema6a, Sema6c, Septin10, Serf2, Sergef, Serhl, Serpine2, Serpinh1, Sesn3, Setd4, Sgcd, Sgsh, Sh3bgrl2, Sh3glb1, Shank1, Shc3, Shc4, Shd, Shf, Shisa6, Shisa8, Siae, Sik2, Sipa1l3, Slc17a6, Slc1a4, Slc24a2, Slc24a3, Slc25a15, Slc25a20, Slc25a5, Slc2a8, Slc35b2, Slc35e2, Slc38a2, Slc38a6, Slc39a9, Slc41a2, Slc41a3, Slc44a1, Slc47a1, Slc66a2, Slc6a11, Slc7a1, Slc8a3, Slc9a1, Slco2a1, Slco2b1, Slco4a1, Slk, Slu7, Smad1, Smarcc2, Smc1a, Smc4, Smc6, Smg7, Smg8, Smurf1, Smurf2, Snhg4, Snph, Snrpd3, Snx32, Socs2, Socs3, Soga1, Sox7, Spaar, Spag4, Specc1l, Src, Srcap, Srek1, Srek1ip1, Srp54a, Srp72, Srrm4, Ssb, Ssbp3, Sst, Ssx2ip, St6gal2, St6galnac2, Stim1, Stim2, Stip1, Stk32b, Stpg1, Strbp, Strn4, Stxbp4, Stxbp5l, Sucla2, Sufu, Sugp2, Sult4a1, Surf2, Swt1, Sympk, Syngr2, Synpo2l, Syt5, Tafa2, Tamm41, Tank, Tbc1d2, Tbc1d22a, Tbc1d30, Tc2n, Tceanc2, Tcerg1, Tecpr1, Tenm2, Terb1, Tesk2, Tet3, Tfam, Tgm2, Thap11, Thoc5, Thoc6, Thop1, Tigar, Tjap1, Tjp1, Tjp2, Tlcd2, Tlcd3a, Tlk1, Tm2d1, Tmc3, Tmcc2, Tmco4, Tmed5, Tmem108, Tmem160, Tmem18, Tmem204, Tmem250-ps, Tmem35a, Tmem35b, Tmem37, Tmem38b, Tmem39a, Tmem39b, Tmem59l, Tmie, Tmpo, Tmtc4, Tnfsfm13, Tnks, Togaram1, Tom1l1, Tomm40l, Top2a, Tor1aip2, Tox2, Tpbgl, Tpd52l1, Tppp3, Trafd1, Trappc2l, Trappc6b, Trim62, Trim71, Trim9, Trip11, Trip12, Trmt2a, Trnp1, Trpc5, Trpm3, Ttc19, Ttc3, Ttpa, Tubd1, Tufm, Tyw3, Uba3, Ubald2, Ube2b, Ube2cbp, Ube2g1, Ube2n, Ube3c, Ubxn2b, Uchl3, Ucp2, Ufc1, Uhrf1bp1l, Ulk3, Umps, Unc119, Unc5d, Urm1, Usp15, Usp3, Usp45, Usp5, Ust, Utp15, Utp25, Uvrag, Vamp8, Vasp, Vcpkmt, Vdac1, Vgll3, Vps13a, Vps13b, Vps26a, Vps4a, Vsig10, Vwa8, Wdhd1, Wdr33, Wdr4, Wdr74, Wdr75, Wdr81, Wfs1, Wipf2, Wnk3, Wrap53, Wrap73, Wtip, Xkr4, Xkr6, Ybx2, Yeats2, Yipf5, Ythdf2, Ywhae, Zbtb1, Zbtb10, Zbtb4, Zbtb46, Zcchc18, Zcchc3, Zfat, Zfhx3, Zfp1, Zfp185, Zfp2, Zfp236, Zfp28, Zfp282, Zfp365, Zfp382, Zfp385b, Zfp407, Zfp414, Zfp472, Zfp516, Zfp523, Zfp609, Zfp688, Zfp703, Zfp768, Zfp770, Zfp777, Zfp786, Zfp787, Zfp853, Zfp867, Zfp94, Zfp971, Zfp994, Zfp995, Zfyve26, Zic3, Zmpste24, Znrd2, Zrsr2, Zswim4,
1110032F04Rik, 1110038B12Rik, 1700029J07Rik, 1700030C10Rik, 1810024B03Rik, 1810059H22Rik, 2610044O15Rik8, 2610306M01Rik, 3110056K07Rik, 4930515G01Rik, 4932441J04Rik, 4933404O12Rik, 4933407K13Rik, 5033430I15Rik, 5330439K02Rik, 6430571L13Rik, 6720427I07Rik, 9930014A18Rik, A230072C01Rik, A230077H06Rik, A430035B10Rik, A830052D11Rik, A930012O16Rik, Abat, Abca13, Abhd5, Ablim3, Abracl, Abt1, Abtb1, Acbd6, Ache, Ackr3, Actr1a, Acvrl1, Adam23, Adamts16, Adamts19, Adamtsl3, Adar, Adarb1, Adck2, Adck5, Adprm, Adra2b, Agap1, Agbl2, Ago3, Agpat2, Aif1, Aig1, Ak4, Akap8l, Aldh18a1, Aldoa, Alkbh8, Aloxe3, Als2cl, Amd1, Anapc1, Anapc15, Angptl4, Ankrd13c, Ankrd33b, Aopep, Ap1s3, Ap3d1, Apba3, Aplp1, App, Araf, Arap2, Arf1, Arf4, Arf5, Arfgef1, Arg2, Arhgap12, Arhgap19, Arhgap30, Arhgap33, Arhgef15, Arid1b, Arid5a, Arl14ep, Arpp21, Arrdc3, Ascc2, Ascl1, Ash2l, Asic1, Asl, Asphd1, Asxl1, Atad2b, Atat1, Atg2a, Atp2b4, Atp5a1, Atp5g3, Atp5pb, Atp6ap1, Atp6v1a, Atp6v1c1, Atp6v1h, Atp9a, Avpi1, B230206L02Rik, B4galt7, B9d1, Babam1, Bahcc1, Baz1b, BC002059, BC005537, BC034090, BC037704, Bcas3, Bcl3, Bdp1, Bend5, Bfar, Bfsp2, Birc6, Blmh, Bloc1s3, Bmp6, Bmt2, Braf, Brca2, Brd3os, Brinp1, Btbd10, Btbd3, Btbd9, Btc, C2, Cacng8, Calcr, Calhm2, Caly, Camk2a, Cant1, Capn1, Card19, Carm1, Casz1, Cav1, Cbr3, Cbx2, Ccdc116, Ccdc117, Ccdc137, Ccdc84, Cchcr1, Cckbr, Ccnd1, Cd276, Cd99l2, Cdc14a, Cdc25a, Cdc25b, Cdk2ap1, Cdk5r1, Cdk5r2, Cdkn1a, Cdon, Cebpd, Cebpzos, Cela1, Celf3, Celf4, Cep135, Cfap36, Cfap53, Cgrrf1, Chchd6, Chdh, Chil1, Chmp2a, Chmp3, Chmp5, Chmp6, Chp1, Chrm3, Chrnb3, Chst14, Chst2, Chsy1, Ciao3, Cic, Cipc, Clcn3, Cldn3, Clip1, Clk2, Clstn1, Cnksr3, Coasy, Cobl, Colec12, Cops3, Cops7b, Cops9, Coq8a, Coro1b, Cox6c, Cplx1, Cracdl, Crb1, Crispld2, Crk, Cry2, Csnk1d, Csnk1g2, Csnk2a1, Cspg4, Ctbp1, Ctnnal1, Ctsa, Cttnbp2, Ctu2, Ctxn1, Cyb5d2, Cyth4, D16Ertd472e, D430040D24Rik, D630045J12Rik, D830050J10Rik, Daam1, Dbndd1, Dctpp1, Dcun1d5, Dcx, Ddi2, Ddn, Ddost, Ddr2, Ddx3x, Ddx5, Ddx6, Degs1, Degs2, Dennd4b, Derl3, Dgcr8, Dgkq, Dgkz, Dhx34, Dhx57, Disp1, Dlec1, Dleu2, Dnaaf9, Dnajc12, Dnajc25, Dnajc6, Dnttip1, Dock8, Dock9, Dot1l, Dpf1, Dusp16, Dusp3, Dync1i2, Dync2i2, E130114P18Rik, Ebf2, Eci1, Ecscr, Edn3, Efcab15, Efhd2, Efs, Egln3, Eid2b, Eif2ak4, Eif2b1, Eif2d, Elf4, Elmo3, Eln, Elovl2, Elovl7, Emc9, Emg1, Endou, Engase, Enoph1, Enpp1, Entpd7, Eogt, Epb41l4b, Epha2, Epha6, Epor, Eps8l1, Etnk1, Evc, Evi5, Exosc8, Eya4, Ezh2, F630040K05Rik, Fam102a, Fam104a, Fam149b, Fam193b, Fam20b, Fam210b, Fam214a, Fam234a, Fam89a, Fasn, Fbl, Fbln7, Fbxl14, Fbxl17, Fbxo21, Fbxo31, Fbxo34, Fbxo41, Fbxo6, Fcgr2b, Fdx1, Fes, Fgf11, Fh1, Ficd, Fkbp11, Fkbp14, Fkbp15, Fktn, Flywch2, Fmnl1, Fmnl3, Fmod, Fmr1, Fndc4, Foxn2, Foxn3, Foxp1, Frat1, Fstl4, Ftsj1, Fut10, Fut8, Fv1, Fxyd6, G3bp2, Gab1, Gabra1, Gabrg3, Gal3st3, Galc, Galnt1, Galnt10, Galnt15, Gapvd1, Garem1, Gas2l3, Gatd1, Gdap1l1, Gdf10, Gfy, Ghitm, Gins2, Gipc1, Gkap1, Glrb, Glrx3, Gm11175, Gm12905, Gm13205, Gm14399, Gm14964, Gm15672, Gm19345, Gm19410, Gm20109, Gm20421, Gm26575, Gm26592, Gm26777, Gm28778, Gm33195, Gm38431, Gm40841, Gm42742, Gm4673, Gm49066, Gm6712, Gm9903, Gm9925, Gmcl1, Gmnn, Gnat1, Gnpda1, Gorasp2, Gosr1, Gpc5, Gpkow, Gprasp1, Grb7, Grik5, Grin1, Grk4, Gsx2, Gtf2f2, Gtf3c6, H1f1, H2-T23, H3f3a, Haus3, Hdgf, Hexim2, Hhex, Hilpda, Hint1, Hltf, Hmgb3, Hmox1, Hnrnpr, Hpf1, Hps4, Hsd17b12, Hsp90b1, Hspg2, Hyi, Ice2, Idh3a, Ifitm1, Ift140, Igf2r, Igfbp4, Igfbp6, Ikzf1, Il10rb, Il15ra, Il17d, Il17rc, Il20ra, Inca1, Ing1, Ing4, Inpp5j, Ints11, Irf4, Ist1, Itga2, Itga8, Itpa, Itpka, Itsn2, Jakmip1, Jrk, Jun, Junb, Katnbl1, Kazald1, Kcna5, Kcnc1, Kcng3, Kcnip1, Kcnj11, Kcnj3, Kcnk1, Kcnma1, Kdm6b, Kl, Klc1, Klf3, Klf4, Klf5, Klhl11, Klhl32, Kmt5a, Kyat1, Large1, Ldlrap1, Letmd1, Lfng, Limch1, Lipt2, Lmod1, Lmtk2, Lpcat4, Lrba, Lrif1, Lrp1b, Lrrc28, Lrrc45, Lrrc7, Lrrc75a, Luc7l, Lypd1, Lypd6b, Lypla2, Macrod1, Maea, Magee1, Mak16, Maml1, Map10, Map3k15, Map3k7, Map6, Marchf7, Marchf8, Mat2a, Mbtd1, Mbtps1, Mbtps2, Mchr1, Mcm6, Mecom, Med12l, Med9, Mef2d, Metap2, Mettl1, Mettl27, Mfrp, Mfsd4a, Mfsd9, Mgat4b, Mgmt, Micu2, Mid1, Mief1, Miga2, Mir8114, Mlxipl, Mmd, Mmp17, Mmp28, Mogs, Morc2b, Mrc2, Mrpl24, Mrpl34, Mrpl58, Mrps30, Mrps31, Msh2, Mta1, Mta3, Mtch2, Mtfp1, Mtg1, Mthfr, Mtif3, Mtor, Mtrfr, Myg1, Myl12a, Myo15, Myo5c, Myo9a, N4bp2l1, Naa30, Naa40, Nadk2, Nampt, Napa, Nars, Nat8l, Nbeal2, Ncam1, Ncaph2, Ncdn, Nceh1, Ncs1, Nde1, Ndor1, Ndufa6, Ndufb2, Necab1, Nek7, Nes, Net1, Nfat5, Nfil3, Nfrkb, Nfya, Nhlh2, Nhlrc1, Nhs, Nkrf, Nln, Nnt, Nop2, Nop56, Npas3, Nphs1, Npr3, Nprl3, Nrde2, Nrip1, Nrn1l, Nsmce1, Nub1, Nudcd1, Nudt18, Nudt22, Nup155, Nup210l, Obi1, Ocrl, Odad3, Odad4, Odf2l, Ogfrl1, Olfml2b, Orc4, Oscp1, Oser1, Otud7b, Otx2, P3h3, Pa2g4, Pafah2, Palb2, Palld, Palm3, Panx1, Papss2, Park7, Parva, Pax2, Pax6, Pcbp3, Pcdh15, Pcdha3, Pcdhga2, Pcgf3, Pcmt1, Pcmtd1, Pde1c, Pde6h, Pde8b, Pdf, Pdia6, Pdk2, Pdlim4, Pdpk1, Per1, Pex11a, Pfkfb4, Pfkp, Pgap4, Pgk1, Pgm3, Phactr1, Phc2, Phf5a, Phkg2, Phospho2, Phyh, Pi4ka, Pi4kb, Pigs, Pik3c2b, Pik3cd, Pik3r4, Pir, Pitpna, Pkig, Pknox2, Plekhb2, Plekhf1, Plekhf2, Plekhm3, Plekho1, Plk2, Plpp2, Pltp, Pmpca, Pnpla2, Podxl, Pola1, Pola2, Poldip3, Polq, Polr1a, Polr3d, Polrmt, Porcn, Ppcs, Ppfibp1, Ppm1e, Ppp1cc, Ppp1r13l, Ppp1r1a, Ppp1r3d, Ppp1r9a, Ppp4r3b, Pramel12, Prcp, Prdm16, Prep, Prkaa1, Prkar2b, Prkcsh, Prkd2, Prkd3, Prox1, Psd2, Psen1, Psma4, Psmc1, Psmd1, Ptbp3, Ptdss2, Ptgs1, Ptms, Ptpa, Ptpn13, Ptpn14, Ptprs, Pum3, Purb, Pxmp4, Qdpr, Rab10os, Rab11fip5, Rab12, Rab13, Rab14, Rab33b, Rab36, Rab3il1, Rabep2, Rad52, Raf1, Rap2c, Raph1, Rasgef1c, Rb1cc1, Rbbp7, Rbbp8, Rbfox1, Rbm25, Rbm26, Rbmx, Rbpms, Rcl1, Rfc3, Rgcc, Rgl1, Rhbdd3, Rheb, Rhno1, Rhof, Ric3, Ripk2, Ripor1, Rmc1, Rnf103, Rnf111, Rnf126, Rnf135, Rnf139, Rnf146, Rnf149, Rnf17, Rnf220, Rnh1, Rnpep, Robo4, Rpap3, Rpl10, Rpl19, Rpl21, Rpl32, Rpl7l1, Rpl8, Rprm, Rps11, Rps26, Rps6ka2, Rps6ka3, Rras, Rreb1, Rrm1, Rsbn1l, Rtraf, Rufy2, Rwdd1, Rxrb, S1pr1, S1pr2, Sbsn, Sc5d, Scamp4, Scara3, Scn4b, Scrn1, Scyl3, Sdf4, Sdhd, Sec62, Sema4g, Senp3, Septin1, Serf1, Serinc2, Serpinb6a, Serping1, Setd1b, Sf3b2, Sf3b4, Sgsh, Sgsm1, Sgsm3, Sh2b1, Sh3bgrl, Sh3bp4, Sh3kbp1, Sh3pxd2b, Sh3rf1, Shc2, Shc3, Shd, Shisa2, Sik2, Simc1, Sin3a, Sin3b, Sipa1, Sipa1l1, Sipa1l3, Slbp, Slc25a1, Slc25a22, Slc27a3, Slc27a4, Slc30a1, Slc35d1, Slc35g1, Slc43a3, Slc4a3, Slc4a8, Slc5a3, Slc7a4, Slc9a3r1, Slc9a5, Slco2a1, Slco5a1, Slit2, Smarcad1, Smim24, Smoc1, Sms, Snrpn, Snx2, Soat1, Socs2, Soga3, Sorcs1, Sos2, Sox12, Sox9, Sp4, Spag1, Spag9, Spata2, Spata2l, Spc25, Spdl1, Spg20, Spg7, Spns2, Spopl, Spryd3, Spsb4, Spx, Sqle, Srbd1, Srd5a3, Srrm3, Srsf11, Srsf7, Ssr4, Sstr1, Sstr4, St6galnac2, St8sia3, Stac, Stambp, Stambpl1, Stat5a, Stat6, Stc2, Stim1, Stradb, Stx17, Sufu, Sult2b1, Sybu, Syf2, Syt4, Tada2a, Tafa1, Tapbp, Tarbp2, Tbc1d1, Tbc1d30, Tbc1d32, Tbc1d5, Tbce, Tbl2, Tbp, Tbrg1, Tcerg1l, Tcp1, Tef, Tfcp2, Thap12, Tiam1, Timm44, Timm8a1, Tipin, Tjp1, Tle1, Tle5, Tlk1, Tlk2, Tll1, Tlr2, Tmcc3, Tmed9, Tmem11, Tmem132e, Tmem145, Tmem150c, Tmem164, Tmem170b, Tmem18, Tmem214, Tmem229b, Tmem39b, Tmem64, Tmem94, Tmem97, Tmsb10, Tnfrsf10b, Tnfrsf12a, Tnfrsf13c, Tnfrsf23, Tnks1bp1, Tnnt2, Tns2, Tob1, Tor1aip1, Tra2b, Trappc2, Trappc6b, Trim39, Trim46, Trim8, Trp53, Trpm1, Trpv2, Trpv4, Tsc2, Tsen2, Tsen34, Tshz2, Tspan11, Tsr1, Ttc13, Ttc21a, Tti2, Ttll11, Ttll3, Ttll7, Ttn, Tub, Tubb3, Tubg2, Tubgcp5, Tufm, Txlng, Txn2, Txndc12, Tysnd1, Uba3, Uba6, Ubap2, Ubc, Ube2f, Ubr1, Ubxn8, Ugp2, Ugt8a, Ulk2, Ulk3, Upf2, Upf3a, Uqcrc1, Urgcp, Usp37, Usp42, Usp48, Utp3, Utp6, Vdac3, Vegfc, Vps54, Vps9d1, Vsir, Vstm2b, Vsx2, Wasl, Wdfy1, Wdr4, Wdr47, Wfs1, Wipf2, Wnk2, Wnk4, Wrap53, Xk, Xkr4, Yeats2, Yipf5, Yod1, Zbtb17, Zbtb21, Zbtb46, Zc3h7a, Zer1, Zfp105, Zfp106, Zfp142, Zfp157, Zfp174, Zfp267, Zfp335os, Zfp354a, Zfp369, Zfp384, Zfp445, Zfp469, Zfp53, Zfp563, Zfp598, Zfp599, Zfp606, Zfp637, Zfp661, Zfp703, Zfp750, Zfp768, Zfp78, Zfp784, Zfp787, Zfp871, Zfp9, Znrf2, Zrsr1,
library(ggvenn)
library(RColorBrewer)
genelist<-subset(merge(tmp, subset(ddply(tmp,.(external_gene_name),nrow),V1>1), by=c("external_gene_name")),external_gene_name!="")
x <- list(
`RTN_HLU+IR`=subset(genelist, group=="RTN_Combined")$external_gene_name,
BRN_HLU=subset(genelist, group=="BRN_HLU")$external_gene_name,
`BRN_HLU+IR`=subset(genelist, group=="BRN_Combined")$external_gene_name,
RTN_HLU=subset(genelist, group=="RTN_HLU")$external_gene_name)
options(repr.plot.width=19, repr.plot.height=20, stroke_alpha=0)
ggvenn(x, show_elements = T, label_sep = ",", fill_color = c("lightblue", "yellow", "lightpink", "grey"),
text_size=6,stroke_color="white") + ggtitle("C")
diffVarRetina<-fread("/home/jupyter/glds202_203/DiffVar/BRN_Combined_vs_Control.tsv",header=T,sep="\t")
subset(corTissue,grepl("Efna5",gene))
subset(diffVarRetina,grepl("Efna5",external_gene_name))
| gene | Type | tissue | Sign |
|---|---|---|---|
| <chr> | <chr> | <chr> | <chr> |
| ensembl_gene_id | ControlvsCombined.mean.log2FC.0 | Statistics.mean | Pvalue.mean | FDR.mean | ControlvsCombined.dispersion.log2FC.0 | Statistics.dispersion | Pvalue.dispersion | FDR.dispersion | external_gene_name |
|---|---|---|---|---|---|---|---|---|---|
| <chr> | <dbl> | <dbl> | <dbl> | <dbl> | <dbl> | <dbl> | <dbl> | <dbl> | <chr> |
| ENSMUSG00000048915 | -0.9195698 | -3.284461 | 0.001021776 | 0.06859507 | -1.749051 | -1.199372 | 0.2303831 | 1 | Efna5 |
tmp<-subset(df, group=="BRN_IR")
tmp$color<-ifelse(tmp$external_gene_name %in% subset(negCorBrain,Type=="IR")$gene, "blue",
ifelse(tmp$external_gene_name %in% subset(posCorBrain,Type=="IR")$gene, "red", "gray"))
mapped = string_db$map(tmp, "external_gene_name",removeUnmappedRows = TRUE )
hits <- mapped$STRING_id
# post payload information to the STRING server
options(repr.plot.width=10, repr.plot.height=10)
payload_id <- string_db$post_payload(mapped$STRING_id, colors=mapped$color )
string_db$plot_network(hits, payload_id=payload_id )
clustersList <- string_db$get_clusters(mapped$STRING_id)
options(repr.plot.width=10, repr.plot.height=10)
par(mfrow=c(1,1))
for(i in seq(5:5)){
tmp<-as.data.frame(clustersList[[i]])
colnames(tmp)<-c("STRING_id")
tmp<-merge(tmp,mapped)
payload_id <- string_db$post_payload(tmp$STRING_id, colors=tmp$color )
string_db$plot_network(tmp$STRING_id, payload_id=payload_id)
}
plot_grid(ggplot(manifest, aes(x=Group,y=Retina_age)) + geom_boxplot(position=position_dodge(0.8),width=0.5) + theme(legend.position="none") + geom_jitter(aes(color=Subject),size=5,position=position_dodge(0.3)) + ggtitle("Retina"),
ggplot(manifest, aes(x=Group,y=Brain_age)) + geom_boxplot(position=position_dodge(0.8),width=0.5) + theme(legend.position="none") + geom_jitter(aes(color=Subject),size=5,position=position_dodge(0.3)) + ggtitle("Brain"))
Warning message: “Removed 1 rows containing non-finite values (stat_boxplot).” Warning message: “Removed 1 rows containing missing values (geom_point).” Warning message: “Removed 3 rows containing non-finite values (stat_boxplot).” Warning message: “Removed 3 rows containing missing values (geom_point).”
library(pheatmap)
corGeneList<-subset(corMatRetina, p.value.Combined<0.05 & estimate.Combined<0 & grepl("Slc",gene))$gene
tmp1<-t(apply(subset(merged203, select=c(subset(manifest,Group=="Combined")$Retina_rna), external_gene_name %in% corGeneList),1,scale))
tmp2<-t(apply(subset(merged203, select=paste0("ratio_",c(subset(manifest,Group=="Combined")$Retina_meth)), external_gene_name %in% corGeneList),1,scale))
tmpMerged<-as.data.frame(cbind(tmp1,tmp2))
colnames(tmpMerged)<-c(paste0(subset(manifest,Group=="Combined")$Subject,"_",subset(manifest,Group=="Combined")$Retina_rna),
paste0(subset(manifest,Group=="Combined")$Subject,"_",subset(manifest,Group=="Combined")$Retina_meth))
labels<-as.data.frame(colnames(tmpMerged))
colnames(labels)<-c("Sample")
labels$Subject<-c(subset(manifest,Group=="Combined")$Subject,subset(manifest,Group=="Combined")$Subject)
labels$Group<-ifelse(grepl("Mmus",labels$Sample),"RNA","RRBS")
labels <- labels %>% column_to_rownames(var="Sample")
tmpMerged<-tmpMerged[ , order(names(tmpMerged))]
rownames(tmpMerged)<-NULL
rownames(tmpMerged)<-paste0(subset(merged203, external_gene_name %in% corGeneList)$external_gene_name,"\n",
subset(merged203, external_gene_name %in% corGeneList)$ensembl_gene_id)
tmpMerged[, "sd"] <- apply(tmpMerged, 1, sd)
tmpMerged<-subset(tmpMerged,sd>0)
tmpMerged<-subset(tmpMerged,select= -c(sd))
dim(tmpMerged)
options(repr.plot.width=6, repr.plot.height=5)
pheatmap(tmpMerged,fontsize = 14,annotation_col=labels,
show_rownames=TRUE,show_colnames=FALSE,cluster_cols=TRUE,
clustering_distance_rows = "correlation",clustering_distance_cols = "correlation",
border_color = NA, main="+ve, Combined")
corGeneList<-subset(corMatMeth, p.value.IR<0.05 & estimate.IR<0)$gene
tmp1<-t(apply(subset(mergedMeth, select=paste0("ratio_",c(subset(manifest,Group=="IR")$Brain_meth)), external_gene_name %in% corGeneList),1,scale))
tmp2<-t(apply(subset(mergedMeth, select=paste0("ratio_",c(subset(manifest,Group=="IR")$Retina_meth)), external_gene_name %in% corGeneList),1,scale))
tmpMerged<-as.data.frame(cbind(tmp1,tmp2))
colnames(tmpMerged)<-c(paste0(subset(manifest,Group=="IR")$Subject,"_",subset(manifest,Group=="IR")$Brain_meth),
paste0(subset(manifest,Group=="IR")$Subject,"_",subset(manifest,Group=="IR")$Retina_meth))
labels<-as.data.frame(colnames(tmpMerged))
colnames(labels)<-c("Sample")
labels$Subject<-c(subset(manifest,Group=="IR")$Subject,subset(manifest,Group=="IR")$Subject)
labels$Group<-c(rep("BRN",length(subset(manifest,Group=="IR")$Brain_meth)), rep("RTN",length(subset(manifest,Group=="IR")$Retina_meth)))
labels <- labels %>% column_to_rownames(var="Sample")
tmpMerged<-tmpMerged[ , order(names(tmpMerged))]
rownames(tmpMerged)<-NULL
rownames(tmpMerged)<-paste0(subset(mergedMeth, external_gene_name %in% corGeneList)$external_gene_name,"_",
subset(mergedMeth, external_gene_name %in% corGeneList)$ensembl_gene_id)
tmpMerged[, "sd"] <- apply(tmpMerged, 1, sd)
tmpMerged<-subset(tmpMerged,sd>0)
tmpMerged<-subset(tmpMerged,select= -c(sd))
dim(tmpMerged)
options(repr.plot.width=5, repr.plot.height=5)
pheatmap(tmpMerged,fontsize = 16,annotation_col=labels,
show_rownames=FALSE,show_colnames=FALSE,cluster_cols=FALSE,
clustering_distance_rows = "correlation",clustering_distance_cols = "correlation",
border_color = NA, main="-ve, IR")
simMatrix <- calculateSimMatrix(as.data.frame(oraGO_Combined)$ID,
orgdb="org.Mm.eg.db",
ont="BP",
method="Rel")
scores <- setNames(-log10(as.data.frame(oraGO_Combined)$`p.adjust`), as.data.frame(oraGO_Combined)$ID)
reducedTerms_Combined <- reduceSimMatrix(simMatrix,
scores,
orgdb="org.Mm.eg.db")
ggplot(reducedTerms_Combined, aes(area = score, subgroup = parentTerm, fill = parentTerm, label = term)) +
geom_treemap(size=0.5, col='white') + ggtitle("Combined") +
geom_treemap_text(col='gray',place="center",reflow=TRUE) +
geom_treemap_subgroup_border(col='white',size=1) +
geom_treemap_subgroup_text(col='white',place="center",fontface="bold",size=12,reflow=TRUE,grow=TRUE,min.size=1) + theme(plot.title = element_text(hjust=0.5,size=10)) +
guides(fill=F)
preparing gene to GO mapping data... preparing IC data... Warning message in calculateSimMatrix(as.data.frame(oraGO_Combined)$ID, orgdb = "org.Mm.eg.db", : “Removed 1 terms that were not found in orgdb for BP” Warning message: “`guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.”